###########################################################
#### figure_1.R                                        ####
#### Generates Figure 1 and accompanying descriptives  ####
###########################################################

#### Calculate state-by-state results ####

res_state = data |> 
  as_survey_design(weights = st_wt_adj) |>
  group_by(state, STATEFIP, state_abb) |> 
  summarise(
    est_affpol = survey_mean(affpol, vartype = "se", na.rm = T),
    est_inparty = survey_mean(inparty, vartype = 'se', na.rm = T),
    est_outparty = survey_mean(outparty, vartype = 'se', na.rm = T),
    est_norm = survey_mean(norm_one_sub, vartype = "ci", na.rm = T),
    est_viol = survey_mean(viol_one, vartype = "ci", na.rm = T)) |> 
  ungroup() |> 
  mutate(fips = str_pad(as.character(STATEFIP), 2, "left", "0"))

#### Generate Maps ####

map_affpol = plot_usmap(data = res_state, regions = 'states', values = "est_affpol") +
  scale_fill_distiller(type = 'seq', palette = 'Purples', direction = 1,
                       name = "Affective Polarization") +
  theme(legend.position = "top")

map_outparty = plot_usmap(data = res_state, regions = 'states', values = "est_outparty") +
  scale_fill_distiller(type = 'seq', palette = 'Reds', direction = 1,
                       name = "Out-Party Affect") +
  theme(legend.position = "top")

#### Generate Distributions ####

mean_df = res_state |> 
  summarise(mean_affpol = mean(est_affpol, na.rm = T),
            cat = "State Mean") |> 
  bind_rows(data |> 
              as_survey_design(weights = wt_adj) |> 
              summarise(mean_affpol = survey_mean(affpol, vartype = 'se', na.rm = T),
                        cat = "Pooled Mean") |> 
              select(mean_affpol, cat))

fig_density = ggplot() +
  geom_density(data = res_state, 
               aes(x = est_affpol, y = ..density.., fill = "State Mean"), 
               alpha = .5) +
  geom_density(data = data, 
               aes(x = affpol, y = ..density.., weight = wt_adj, fill = "Pooled Mean"), 
               alpha = .5) +
  geom_vline(data = mean_df, aes(xintercept = mean_affpol, color = cat), show.legend = F) +
  scale_fill_manual(name = "Data",
                    values = c("#440154FF", "#22A884FF"),
                    labels = c("Pooled", "State")) +
  scale_color_manual(values = c("#440154FF", "#22A884FF")) +
  theme_prl() +
  labs(x = "Affective Polarization", y = "Density") +
  theme(legend.position = c(.2,1.2))

#### Combine Plots and Export ####

fig_1 = map_affpol + map_outparty + fig_density

print(fig_1)

ggsave(here("Plots", "figure_1.pdf"), fig_1,
       units = "in", height = 4, width = 10, dpi = 600)

#### Descriptives ####

## Min/Max ##
print("Maximum, Minimum, and Quartiles for Outparty Affect")
slice_max(res_state, order_by = est_outparty) |> select(state, est_outparty, est_outparty_se)
slice_min(res_state, order_by = est_outparty) |> select(state, est_outparty, est_outparty_se)
quantile(res_state$est_outparty, c(.25,.5,.75))

print("Maximum, Minimum, and Quartiles for Affective Polarization")
slice_max(res_state, order_by = est_affpol) |> select(state, est_affpol, est_affpol_se)
slice_min(res_state, order_by = est_affpol) |> select(state, est_affpol, est_affpol_se)
quantile(res_state$est_affpol, c(.25,.5,.75))

## % of Total Range ##
print("% of Total Range")
(59.9-42.1)/201 

## ICC ##
between_var = data |> 
  group_by(state_abb) |> 
  summarise(mean_affpol = mean(affpol, na.rm = T)) |> 
  summarise(bvar = var(mean_affpol, na.rm = T)) |> 
  pull(bvar)

within_var = data|> 
  group_by(state_abb) |> 
  mutate(mean_affpol = mean(affpol, na.rm = T),
         centered_affpol = affpol - mean_affpol) |> 
  ungroup() |> 
  summarise(wvar = var(centered_affpol, na.rm = T)) |> 
  pull(wvar)

state_icc = between_var / (between_var + within_var)
print("ICC")
state_icc

#### Clean Up ####
rm(res_state, map_affpol, map_outparty, fig_density, between_var, within_var, state_icc)
